library(tidyverse)

── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr     1.1.4     ✔ readr     2.1.5
✔ forcats   1.0.0     ✔ stringr   1.5.1
✔ ggplot2   3.5.1     ✔ tibble    3.2.1
✔ lubridate 1.9.4     ✔ tidyr     1.3.1
✔ purrr     1.0.2     
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag()    masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(lubridate)
library(tsibble)

Registered S3 method overwritten by 'tsibble':
  method               from 
  as_tibble.grouped_df dplyr

Attaching package: 'tsibble'

The following object is masked from 'package:lubridate':

    interval

The following objects are masked from 'package:base':

    intersect, setdiff, union

library(feasts)

Loading required package: fabletools

library(dataRetrieval)
library(plotly)

Attaching package: 'plotly'

The following object is masked from 'package:ggplot2':

    last_plot

The following object is masked from 'package:stats':

    filter

The following object is masked from 'package:graphics':

    layout

poudre_flow <- readNWISdv(
  siteNumber = "06752260",
  parameterCd = "00060",
  startDate = "2013-01-01",
  endDate = "2023-12-31"
) |>
  renameNWISColumns() |>
  mutate(Date = yearmonth(Date)) |>
  group_by(Date) |>
  summarise(Flow = mean(Flow, na.rm = TRUE)) |>
  ungroup()

GET:https://waterservices.usgs.gov/nwis/dv/?site=06752260&format=waterml%2C1.1&ParameterCd=00060&StatCd=00003&startDT=2013-01-01&endDT=2023-12-31

# Convert to tsibble
poudre_ts <- as_tsibble(poudre_flow, index = Date)

# Base plot
flow_plot <- ggplot(poudre_ts, aes(x = Date, y = Flow)) +
  geom_line(color = "steelblue") +
  labs(title = "Monthly Average Streamflow (Cache la Poudre River)",
       x = "Date",
       y = "Discharge (cfs)") +
  theme_minimal()

# Animate with plotly
ggplotly(flow_plot)

# Subseries plot
poudre_ts |>
  gg_subseries(Flow) +
  labs(title = "Seasonal Subseries Plot of Streamflow",
       y = "Discharge (cfs)")

In the subseries plot, each line represents the monthly averages across all years for that month.Seasons are the 12 calendar months, and each point shows the average flow in that specific month over time. The subseries help us visualize seasonal variation—peaks are typically in late spring/early summer due to snowmelt, while lows are in winter.

# STL decomposition
stl_model <- model(poudre_ts, STL(Flow ~ season(window = "periodic")))

# Visualize decomposition
components(stl_model) |>
  autoplot() +
  labs(title = "STL Decomposition of Streamflow Data")

The STL plot shows three parts:

Trend: The long-term change in streamflow. There’s some fluctuation over the years, maybe tied to precipitation variability.

Seasonal: Recurring annual patterns, peaking in spring/summer.

Remainder: Unexplained variability (like extreme weather events or human activities like dam releases).

The seasonal component is very strong, showing snowmelt influence. The trend fluctuates mildly showing some drought impacts or climate variability.